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We investigate the oscillations of slowly rotating superfluid stars, taking into account the vortex mediated mutual friction 
force that is expected to be the main damping mechanism in mature neutron star cores. Working to linear order in the 
rotation of the star, we consider both the fundamental f-modes and the inertial r-modes. In the case of the (polar) f-modes, 
we work out an analytic approximation of the mode which allows us to write down a closed expression for the mutual friction 

damping timescale. The analytic result is in good agreement with previous numerical results obtained using an energy 
integral argument. We extend previous work by considering the full range of permissible values for the vortex drag, e.g. the 
friction between each individual vortex and the electron fluid. This leads to the first ever results for the f-mode in the strong 
drag regime. Our estimates provide useful insight into the dependence on, and relevance of, various equation of state 
parameters. In the case of the (axial) r-modes, we confirm the existence of two classes of modes. However, we demonstrate 
that only one of these sets remains purely axial in more realistic neutron star models. Our analysis lays the foundation for 
companion studies of the mutual friction damping of the r-modes at second order in the slow-rotation approximation, the first 
time evolutions for superfluid neutron star perturbations and also the first detailed attempt at studying the dynamics of 
superfluid neutron stars with both a relative rotation between the components and mutual friction. 

I. INTRODUCTION 

Neutron stars have a complex interior structure. With core densities reaching several times the nuclear saturation 
density, these objects require an understanding of physics that cannot be gained from laboratory experiments. This 
makes the modelling of neutron star dynamics an interesting challenge. On the one hand, one has to consider 
exotic physics that is, at best, poorly constrained. On the other hand, one may ask to what extent observations 
can distinguish between different possible models. An excellent example of this interplay concerns the possibility 
that the quarks may deconfinc in the high density region. If this is the case, it will have a considerable effect on 
transport properties associated with viscosity and heat conductivity. In fact, such a quark core is expected to be 
a colour superconductor [l|. The dynamics of this exotic state of matter, and the relevance of its different possible 
phases, is not yet certain. In order to improve our understanding of this problem, we need to build more precise 
stellar models and study, for example, their oscillation properties in detail. In this context, considerable attention 
has been focused on the inertial r-modes of a rotating star. The r-modes are interesting because they can be driven 
unstable by the emission of gravitational radiation, see 0, 0] for literature reviews. The r-modc instability window is, 
however, sensitive to the physics of the neutron star interior. Since the bulk and shear viscosities are quite different 
in a quark core, compared to "normal" npe matter, one may hope to use observations to constrain the theory, see 
[l[ for a discussion of the relevant literature. In absence of a direct detection of an r-mode gravitational-wave signal, 
this analysis would have to be based on the nature of the instability window. The idea would be that, if an observed 
neutron star spins at a rate that would place it inside a predicted instability region, one may be able to rule out this 
particular theoretical model. Of course, this argument comes with a number of caveats. It could, for example, be 
that additional physics places a stronger constraint on the r-modes than the considered mechanisms. Inevitably, this 
becomes a "work in progress" where improved theoretical models are tested against better observational data. 

In order to consider "realistic" neutron stars, it is important to appreciate the relevance of superfluidity. A neutron 
star is expected to contain a number of superfluid/superconducting components and it is crucial to understand 
to what extent this affects the stars oscillation properties. It is well established that the behaviour of a superfluid 
system can differ significantly from standard hydrodynamics. The most familiar low-temperature system is, perhaps, 
He4, which exhibits superfluidity below a critical temperature near 2 K. Experimentally, it has been demonstrated 
that this system is very well described by the Navicr-Stokcs equations above the critical temperature. Below the 
critical temperature the behaviour is very different, and a "two-fluid" model is generally required (see Q for a very 
recent discussion). Superfluid neutron stars are, to some extent, similar. The second sound in Helium is analogous to 
a set of, more or less distinct, "superfluid" oscillation modes @, B ^ neutron star. These additional modes arise 
because the different components of a superfluid system are allowed to move "through" each other. The dissipation 
channels in a superfluid star are also quite different. Basically, the superfluid flows without friction. In the outer core 
of a neutron star, which is dominated by npe matter, one expects the neutrons to be superfluid while the protons 
form a superconductor. As a result, the shear viscosity is dominated by e-e scattering The bulk viscosity, 

which is due to the fluid motion driving the system away from chemical equilibrium and the resultant energy loss due 
to nuclear reactions, is also expected to be (exponentially) suppressed in a superfluid [3]. These effects have direct 
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implications for the damping of neutron star oscillations, and play a key role in determining the r-mode instability 
window for a mature neutron star. This is, however, not the end of the story. A superfluid exhibits an additional 
dissipation mechanism, usually referred to as "mutual friction" . The mutual friction is due the presence of vortices 
in a rotating superfluid. In a neutron star core, the electrons can scatter dissipatively off of the (local) magnetic field 
of each vortex (see [HI, for discussions and references). This effect may dominate the damping of realistic 
neutron star oscillation modes. 

The basic requirements of a rudimentary model for superfluid neutron star oscillations should now be clear. One 
must account for the additional dynamical degree(s) of freedom, and also account for the mutual friction damping. 
This is obviously only a starting point, but the problem is sufficiently complicated that one may want to proceed with 
care. There has already been a number of studies of dissipative superfluid oscillations. The area was pioneered by 
Lindblom and Mendell who, in particular, demonstrated that the gravitational-wave instability of the fundamental 



f- modes would be suppressed in a superfluid star Following the discovery of the r-mode instability, they also 



provided the first accurate estimates of the relevance of the mutual friction for these modes [l5l | . Similar results were 
subsequently obtained by Lee and Yoshida • These studies provide important assessments of the relevance of the 
mutual friction damping. There are, however, a number of reasons why we need to return to this problem. Most 
importantly, we want to consider more realistic neutron star models, including finite temperature effects, magnetic 
fields and the possible presence of exotic (hypcron and/or quark) cores. The additional physics brings additional 
complications, like additional fluid degrees of freedom, boundary layers at phase-transition interfaces and fundamental 
issues concerning dissipative multifluid systems [T3| . We also need to move away from the assumption that the 
vortex drag, which leads to the mutual friction, is weak. Strong arguments suggest that this is not going to be 
the case when the protons form a type II superconductor and there are magnetic fluxtubes present in the system 
0, IHj HOl- The neutron vortices may be "pinned" to the fluxtubes leading to a strong drag regime. The strong 
drag problem has only been considered recently (2ll . [2^ . and the first results demonstrate the presence of a new 
instability in systems where the two cornponents rotate at different rates. This instability, which may be relevant 
for the understanding of pulsar glitches [24I, provides a direct demonstration that the dynamics in the strong drag 
regime may be both complicated and interesting. The present investigation lays the foundation for future work in 
this direction by allowing for strong drag. In particular, we retain the dynamic contribution to the mutual "friction" 
that has previously been neglected as a matter of course. 



II. THE TWO-FLUID EQUATIONS 

Our discussion is based on the standard two-fluid model for neutron star cores [13, H^- That is, we consider two 
dynamical degrees of freedom loosely speaking representing the superfluid neutrons (labeled n) and a charge-neutral 
conglomerate of protons and electrons (labeled p). Assuming that the individual species are conserved, we have the 
usual conservation laws for the mass densities px, 

atPx + V,;(p,<) =0 (1) 

where the constituent index x may be either p or n. Meanwhile, the equations of momentum balance can be written 

{dt + viV, ) [vf + s.wr) + V, (Ax + <&) + ex< V,«!^ = /f /px (2) 

where the velocities are v^, the relative velocity is defined as lu^y = v"^ — Vy and p.^ — Mx/jtIx represents the chemical 
potential (we will assume that mp = mn throughout this paper). $ represents the gravitational potential, and the 
parameter ex encodes the non-dissipative entrainment coupling between the fiuids [l7l.[23t. The force on the right-hand 
side of ^ can be used to represent various other interactions, including dissipative terms fv^ - 

In the following we will focus on the vortex-mediated mutual friction. Assuming that the two fiuids exhibit solid 
body rotation we have a force of form [l^ (see also [ll|, [T^j ) 

/f = 2pnB'e,jkn^w^y + 2p^Be,,kn'e'''"'niw^ (3) 

Here, fl^ is the angular frequency of the neutron fluid (a hat represents a unit vector). The mutual friction parameters 
are intimately related to the induced friction on the vortex. The latter is often described in terms of a dimensionless 
"drag" parameter TZ such that 

B' = nB= — — ^ (4) 

In the standard picture, the mutual friction is due to the scattering of electrons off of the array of neutron vortices. 
This leads to 7?. <C 1, i.e., B' ^ B, and hence the first term in the mutual friction force can be ignored. There are. 
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however, good arguments for why the problem may be in the opposite regime. In particular if one considers the 
interaction between the fluxtubes in a type II proton superconductor and the neutron vortices [l^, [l^ . [20| . Then one 
would expect to be in the strong drag regime where 7Z ^ 1, i.e., B' ^ 1 while B remains small. Superfluid oscillations 
in this regime have not previously (with the exception of [22| ') been considered. 

Anyway, from ([3]) wc sec that the mutual friction will not be present in a non-rotating star. This is obvious since 
there would then be no vortices in the first place. Of course, any non-trivial motion of the superfluid neutrons leads 
to vortex generation. This means that a generic perturbation of a non-rotating star will be associated with a local 
vorticity which could lead to mutual friction. However, in this context the resulting mutual friction interaction would 
require a perturbative calculation to be carried out to second order. As far as we are aware, such calculations have 
not yet been attempted. It may be an interesting problem for the future. 



III. THE PERTURBATION EQUATIONS 

A. Decoupling the degrees of freedom 

If we want to consider the effects of mutual friction we need to consider rotating stars. To keep the problem 
tractable (at least initially) we assume that the background configuration is such that the two fluids rotate together. 
Perturbing the equations of motion and working in a frame rotating with Cl^ we then have 

dt{6v^ + ejwj'') + W,{Sfi. + 5$) + 2e,,kn^Sv^^ = S{f^/p^) (5) 

and 

dtSp^ + V,{p^Svi) = (6) 

where 6 represents an Eulerian variation. 

From previous work on superfluid neutron star perturbations (and indeed the large body of work on superfluid 
Helium) we know that the problem has two " natural" degrees of freedom, see for example [y, 0, H, [13, [1^, [1^ . One 
of the degrees of freedom represents the total mass flux. Introducing 

pSv^ = PnSv^^ + PpSv^p (7) 

and combining the two Euler equations wc find that 

dtSv^ + Vj(5$ + -\7,5p - ^5p\/ip + 2e,jk^^5v'' ^ (8) 
P P 

where p = Pn + Pp and the pressure is obtained from [sj 

ViP = Pn^tPn + PpVi/ip (9) 

In deriving this relation we have used 

PnVi(5/2n + Pp^iSpp = Vi(5p - 5pS/ijl = Vi(5p - ^5p\7^p (10) 

where it has been assumed that the two fluids are in chemical equilibrium in the background. That is, we have 
pn = Pp = P- We also have the usual continuity equation 

dtSp + Vj{pSv^) =0 (11) 

At this point we have two equations which are identical to the perturbation equations for a single fluid system. It 
is particularly notable that ^ does not have a force term. This follows immediately from the fact that we are only 
considering the mutual friction interaction. In other situations, say including shear viscosity, we would no longer have 
a homogeneous equation. 

Of course, we are considering a two-fluid problem and there is a second degree of freedom to take into account. To 
describe this, it is natural to consider the difference in velocity. Thus, we introduce 



Sw-' — dv^ — Sv^ 



(12) 
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Combining the two Euler equations in the relevant way we have 

(1 - e)^t6w^ + ViSp + 2B'eijk^Hw^ - Be.jk^' e^^'^niSw^n = (13) 

Here we have defined 

5(3 = Sjlp - 5fLn (14) 

which represents the (local) deviation from chemical equilibrium induced by the perturbations. We have also intro- 
duced the simplifying notation 

s = en/xp, B' ^l~B'/xp , B^B/xp (15) 

where Xp ^ pp/ p is the proton fraction. Again, equation (|13[) does not couple the different degrees of freedom. 

The coupling is entirely due to the second continuity equation. It is natural to use the proton fraction to complement 
the total density p. Then we find that 

dtSxp + -Vj [a;p(l — Xp)pSw^] + Sv^V jXp ~ (16) 

This equation shows that the two dynamical degrees of freedom are ex plic itly coupled unless the proton fraction is 



constant. This fact has already been pointed out by Prix and Rieutord [24 1 



Before moving on, it is worth asking to what extent it is possible to find solutions that are purely co-moving, i.e. 
for which 5w^ =6/3 = 0. From the above equations it is easy to see that such a solution would have to satisfy 

dtSxp + Sv^VjXp = (17) 

This condition is trivially satisfied if the proton fraction is uniform. In addition, it will be satisfied for fluid motion 
that has (for a spherical background configuration) no radial component and also do not lead to variations in 6xp. 
Are there oscillation modes with this character? Indeed, to leading order in the slow-rotation approximation the 
standard r-mode satisfies these criteria. It is purely axial and the associated density perturbations appear at order 
Q,^. However, in general we do not expect to find any oscillations of a "realistic" neutron star model to be purely 
co-moving. This means that a generic neutron star oscillation mode will be affected by mutual friction. 



B. Boundary conditions 

To completely specify the perturbation problem, we need boundary conditions. At the centre of the star we simply 
require that all variables are regular. The surface of the star is somewhat more complex. In reality one does not 
expect the superfluid region to extend all the way to the surface. A real neutron star will always be covered by a 
single fluid envelope (eg. the outer parts of the elastic crust). However, for simplicity we do not want to deal with the 
various interfaces in the present analysis (isl. [27l. [28j . Instead we will consider stars with a two-fluid surface, which is 
obviously somewhat artiflcial. 

A reasonable approach is to assume that the perturbed star has a unique surface. That is, let the two perturbed 
fluids move together (in the radial direction) at the surface. In the two-fluid problem we have two distinct Lagrangian 
displacements [2^. These follow from 

dtc. = svi+^iv,vi+viv,c. (18) 

We have assumed that the two fluids corotate in the background, i.e. we have v^^ = Vp. If we also impose that there 
is a common surface, then we have ^"^^ = f p at r = i? and it follows that we should require 

Sw'' = Svl, - (5< , at r = i? (19) 

From (jl3p we see that this implies that, for a non- rotating configuration we must also have 

drSp = , at r = i? (20) 

When we determine the rotational corrections to the f-mode we will impose this condition also at first slow-rotation 
order. This is not entirely consistent, cf. but it is straightforward to relax this condition should it be required. 

If the two fluids move together at the surface, it also follows that 

Sp + p^^VjP = Sp + C'drP = Ap = 0, atr = R (21) 

where = = ^p at the surface. This is the usual single fluid condition of a vanishing Lagrangian pressure variation 
Ap. 
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C. A bit of chemistry 

Consider the various equations that we have written down. At this point the two degrees of freedom [Svi , 6p\ and 
[Swi^ SP] only couple explicitly through p6p . In fact, if we assume that the two fluids are incompressible, then there is 
no coupling at all. Since the mutual friction only enters the problem via (jl3p . it is thus the case that any incompressible 
dynamics in the [6vi,5p] sector will be unaffected by mutual friction. This shows that, if we are interested in the 
effect of mutual friction on (say) the f-mode oscillations of a star it is not meaningful to consider an incompressible 
model. We know already from the outset that we would only find the mutual friction effects on the counter-moving 
"superfluid" modes. This may be an interesting problem, but it is not our main motivation here. 

For compressible models, the two degrees of freedom also couple indirectly. Basically, we need to use the equation 
of state to relate [Sp, 6/3] to [Sp, Sxp] . For models where the two fluids co-rotate in the background we can use l3Sj 

and 

Using these relations, or their "inverse" , we see that the two degrees of freedom couple in a more subtle way. If we 
choose to reduce the problem by eliminating Sp and 6f3 then the coupling arises through the boundary conditions and 
the Euler equations. If, on the other hand, we eliminate 6p and Sxp then the coupling enters through the continuity 
equations. 



IV. DISSIPATION INTEGRALS 



In order to estimate the damping associated with various dissipation mechanisms one can either (i) account for the 
dissipative terms in the equations of motion and solve for the damped modes directly, or (ii) solve the non-dissipative 
problem and use an energy integral argument to estimate the damping rate. In the typical situation when the damping 
is very slow the second strategy should be reliable. Indeed, all studies of damped neutron star oscillations have used 
this approach (see Q for a discussion). Given this, it is natural to pause and consider the energy integral approach 
to the mutual friction problem. 



A. The conserved energy 

To work out a suitable energy associated with a given perturbation, we first multiply (O with Px<5wx (where the 
bar represents complex conjugation). Then we add the result to its complex conjugate. Combining the individual 
contributions from the neutron and proton fluids and integrating over the star we find that, when /f = 0, the result 
is a total time derivative of two terms. The first term is the "kinetic energy", which follows from 

Ek = \j [{pn - 2«)|(Si'„P + AaRc{5vl5v^^) + (pp - 2a)\6vp\''] dV (24) 

where 2a = Px^x- Alternatively, expressing this in the co- and countermoving variables, we have 

^k = ljp [M'' + (1 - £)2^p(l - ^p)\H^] dV (25) 

The "potential" energy requires a bit more work. Using the divergence theorem and the continuity equation one 
can show that we need 

Px[KV»((5px + <^*) + c.c]dV = J pSvliSp^ + S^) + c.c]eldS + j [{Sjly, + 5^)dt5p^ + c.c]dV (26) 

(c.c. represents the complex conjugate). The surface term vanishes if px ^ at the surface of the star. It also 
vanishes for modes that have no radial component. Adding the contributions for the neutron and proton fluids we 
see that we need 



/ 



5^dt5pdV = -J— [ S^dtV^S^dV = / 6^el.V,dt6^dS - [ (V,<5$)(9t(V' dtS^)dV (27) 
47rG J 4ttG J AttG J 
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Again, one can argue that the surface term vanishes. 
Finally, we have terms of form 

Sfl^dtSp^dV 

It is natural to express these in terms of the perturbed densities using (this is valid for co-rotating background models 
only) 



Adding all the terms together we find that the "potential energy" follows from 



Er, 



1 



dpn 



Pp 



dpp 



Rc (SppSpn) 



dPvJ p. 



|(V,<5$)| = 



dV 



(28) 



(29) 



or 



dp 

8(3 



[2Re{5h5(3) 



With these definitions it follows that the total "energy" is conserved, i.e. 

dtE = dt{Ek + Ep) ^ 

when /f = 0. These energy expressions are equivalent to those used by Lindblom and Mendell ^ 



(30) 



(31) 



B. Mutual friction 



Even though we will include the mutual friction terms in the equations of motion, it is useful to work out the 
correspon ding dissipation integrals. After all, this is the way that the mutual friction damping has traditionally been 
evaluated |l4l. [isl. [l6l| and we want to be able to compare the two approaches. 

First consider the B' terms. It is easy to show that these terms are not dissipative. We find that 



2dtEB' = 2 / B'e,,k^^ 



(5<(5w„p + (5up(5wp„ + c.c. 



dV ^0 



(32) 



by symmetry. This result is not surprising. In fact, we see from (|13p that the B' terms enter the equations of motion in 
the same way as the Coriolis force. Since the Coriolis terms vanish identically when we multiply each Eulcr equation 
with this should be true also for the non-dissipative part of the mutual friction. 
Finally, it is straightforward to show that the dissipative terms lead to 



dtEts 



PnBiSlS'J" - SrS'j]Cimi[SvlSw^P + (55;5<^ + c.c.]dV 



p,m[5T - n'''n^\5wlJwZdV (33) 



Let us now ask how we can use these results to estimate the mutual friction damping timescale. Let us assume that 
we have a mode solution to the full dissipative problem. That is, we have a solution with time dependence e*'^* where 
Lj = Ljj. -\- i/r such that r is the damping timescale. From the fact that the energy is quadratic in the perturbations 
it follows that |35| 



2E 



dtE 

Moreover, since the solution satisfies the dissipative equations of motion we also know that 

dtE = dtEe 

Hence, we can equally well use 

2E 



(34) 



(35) 



(36) 
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As long as we are using the complete solution to evaluate this expression, it is an identity. However, in many cases 
we do not have access to the solution to the dissipative problem. (If we did, we would not need the energy integrals 
in the first place.) In these cases we can still estimate the damping timescale by evaluating the right-hand side of 
pop using the non-dissipative mode solution. When the damping is sufficiently slow, in the sense that the dissipative 
terms have a small effect on the eigenfunctions, this estimate should be reliable. Of course, one should not expect it 
to yield exactly the same result as the solution to the full dissipative problem. 



C. Gravitational- wave emission 



Finally, let us work out the multipole formulas for gravitational-wave emission from a two-fluid star. This exercise 
is particularly relevant if we are interested in oscillations that may be driven unstable by gravitational-wave emission 
. The main motivation for including it here is that it demonstrates the intuitive result that gravitational waves are 
only generated by the co-moving degree of freedom. 

Following [30| we need the mass multipoles 



5Dim = / TooYi^r'dV « / STo^Yi^r'dV (37) 

and the current multipoles 

5Jim - / {-ro,)Y,%,dV « / {~5T^^)Y,%dV (38) 



In these expressions Yim are the standard spherical harmonics and Y^i^ are the magnetic multipoles |3(3| . 

To work these out we start with the usual expression for the two-fluid stress-energy tensor in general relativity [111 

T^. = *.9p. + + nPfiP (39) 

In the relativistic formulation, see (sij for a review and a survey of the literature, the central variables are the fluxes 
n^. The associated momenta follows from 

^il^B^nl + A-ynl (40) 

Hence, the A^^ coefficients encode the entrainment effect. Let us now work in the frame of an observer moving with 
four-velocity such that 

< = [7x, 7x<] , with 7x - (1 - vl)-'^^ (41) 
where v^. is the associated three- velocity. Then it follows that 

Too = ^-500 + nlB'^jl + 2n„np^"P7„7p + nlBP-/l (42) 
We want the Newtonian (low-velocity) limit of this expression. Thus we let 

500 -1 , and 7x 1 (43) 

and we get 

Too « + nlB^ + 2nnJip^"P + ji^^p (44) 

Finally use the definition 

* = A - ny^^ - nlnl « A + n^S" + 2n„npyl"P + n^gP (45) 

and A = —p to arrive at 

Too= P — > Too ~ 5p (46) 

as one would have expected. 

For the current multipoles wc need 

To, = nlB'^jlvJ + n„7.p^"P7„7p(«J + vj) + nlB^^^v^ (47) 
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In the low- velocity limit, this leads to 

To, » n„(nnS" + np^"P)Wj" + np{7ipBPv] + n„^"P)t-P (48) 

Finally, we need to write this expression in terms of the Newtonian variables. This can be done by comparing the 
momenta, 

^ « + ^^vl = (1 - £n)< + e^v^ (49) 

m m m ^ ^ 

This suggests that we should identify 

= ! — £„, and = £„ (50) 

m m 

Using analogous expressions for the protons we see that 

Toj ~ PnVf + PpV^ (51) 

which leads to (for a co-rotating background) 

Toj ~ PnSvj + PpSv^ + Spvj = pSvj + Spvj (52) 
These results show that it is only the co-moving degree of freedom that radiates gravitationally. 

V. SLOW ROTATION PERTURBATION EQUATIONS 

Let us now return to the problem of oscillating superfluid neutron stars. We will first derive the general perturbation 
equations for a slowly rotating superfluid star. To do this we expand all variables in spherical harmonics. Since we 
expect rotation to couple the various multipoles, we represent the velocity perturbations by the general expressions 



and 



1 /I Til \ 1 

r V sm d I r'^ sm d ^ 



(53) 



(54) 



Note that we represent the "co-moving" degree of freedom by the uppercase amplitudes [M^;, V;, Ui\ while the "counter- 
moving" degree of freedom corresponds to the lowercase quantities All scalar perturbations are expanded 
in spherical harmonics, i.e. we have 5p = Spif^ etcetera. From now on the sum over I will be implied. 

One can use a number of different strategies in writing down the perturbation equations. To some extent this is a 
matter of taste. However, in the slow-rotation problem it can be advantageous to work with a set of equations where 
the coupling between different multipoles is minimal. The set of equations that we use was chosen using this criterion. 
We also decided to use the velocity perturbations as our main variables. This approach is analogous to that used by 
Lockitch and Friedman [s^l in their analysis of inertial modes of single fluid stars. It is notably different from the 
two-potential formalism pioneered by Ipser and Lindblom [ssj, which was extended to superfluid stars by Lindblom 
and MendeU [H. 

We replace each of the perturbed Euler equations with three equations. The first is the radial component of the 
vorticity equation that follows if we take the curl of or (|13p . Assuming that the perturbations have a harmonic 
dependence on time, exp(ici;t), we get 

+ - 2mr!]i7,r/" + + 2){Wi - Wi)Qi+xYl'^^ - 2^(1 - l)[Wi + {I + l)l^]Q/y,"_\ = (55) 

and 

{l{l + l)cj(l - e) - 2mnB' - 2in[l{l + 1) - rr?]B} uiYi"" 
+2Q[{1 + 2)B' - imB] {wi - Ivi} Qi+iYi^^^ 
-2n[{l-l)B' + i7nB]{{l + l)vi+wi}QiY/I^^^0 (56) 
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In deriving these equation we have made use of the standard recurrence relations 

cos err = Qi+iYi^i + QiYi^i (57) 

and 

sin edeVr = lQi+iYi"^i - (l + l)QiYi^i (58) 

where 



Qi 



{21- 1)(2/ + 1) 



(59) 



For future reference, note that Q„i = and Q^+i = l/(2m + 3). 

Next we could use also the 9 (or ip) components of the vorticity equation. However, as discussed in [s^l there is 
a slightly simpler alternative. We first of all use a pair of equations analogous to the "divergence" equation in [33 |. 
These can be written 

{[l{l + l)uj - 2mn]Vi - 2mnWi - + 1)[(1 - Xp)5p^ + a;p<5/if]} Yj" 

-2m{l + 2)C/zOz+ir/;\ - 2n{f - 1)UiQiYi"1i = (60) 

and 

- + 1)5/3/ + 2mQ.B'wi - + l)w(l - e) - 2mnB']vi 

+2inB + i)Q'i+i + 2)Qf] wi 

+2iW [m^ + l{l + 3)gf+i + 2)Qf] vi} Yi"' 

-2n{l + 2){IB' + imB)uiQi+iYi'^^ - 2^l{l ~ + 1)6' - Mn6]u/Q,Y,"i 

-2inB{l + i){lvi ~ wi)Qi+iQi+2Yi"^^ ~ 2inB{l - 2){wi + (^ + \)vi\QiQi^^Y{'l^ = (61) 

Meanwhile the radial components of the Euler equations lead to 

- Xp)rdrSilf + XprdrSfl^i] + 2mnVi - uWijYi"' 
+ 2mUiQi+iYi'!^^ - 2n{l + l)UiQiYi'!^^ = (62) 

and 

{irdrSpi - [c^(l - e) - 2in{l -Qf- Qf+i)Bwi + 2n[mB' - + l)Qj - lQ^i+^)B]vi] 
+2n{lB' + imB)uiQi+iY{^^ - 2n[{l + 1)6' - imB]uiQiYi'^^ 

+2in{lvi - wi)BQi+iQi+2Yi'^^ - 2ir![(? + l)vi + wi]BQiQi^iY{^2 = (63) 
Finally, the continuity equations become 

iujr^5pi+dr{rpWi)-l{l + l)pVi^Q (64) 

and 

iuipr^Sxi + dr[xp{l — Xp)rpwi] — Xp(l — Xp)l{l + l)pvi + pWirdrXp = (65) 

This completes the description of the general first order slow-rotation problem. 

In order to deduce the relevant recurrence relations from the above equations we need to recall that we have been 
implying summation over I. That is, we are considering relations of form 

^ [aiQiQi-iY/!^^ + biQiY/!^, + ciYr + diQi+iYi"^, + eiQi+iQ,+2Yl^^] = (66) 
I 

Using orthogonality of the spherical harmonics, i.e. multiplying by Y"^ and integrating over the sphere, we obtain 
the recurrence relation 

Qn-lQn = (67) 

Given this result, it is straightforward to write down recurrence relations for the various classes of oscillation modes 
of a rotating superfluid star. However, since the level of rotational coupling is different for different kinds of modes, 
it is not particularly useful to write down the general relations. Instead, we focus on two specific examples. 
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VI. THE F-MODES 

Let us begin by considering modes that are non-trivial already in a non-rotating star. Then we first need to solve 
the non-rotating (and non-dissipative since the mutual friction damping requires rotation) problem. Simply setting 
n = in our perturbation equations we see that the polar and axial degrees of freedom decouple (as they should) . It 
is also clear, cf. (|55|) and ((56)) . that there will not exist any purely axial modes in the non-rotating case. This means 
that we can make the Ansatz 

U! ^ LOq + LOlQ (68) 

together with 

Wi = + nwl , Vi^Vi° + nVi^ , Ui ^ wl (69) 

and 

wi = wi + flwl , vi = v'^ + rivl , ui = nu\ (70) 

and similarly for the various scalar perturbation quantities. For example, in the case of the proton fraction wc have 
fep = Y.I ^xiYl^ with 

5x1 = 5x°i + nSxj (71) 

A. The non-rotating problem 

At the non-rotating level the equations in Section IV provide the following relations 

PnSjlli + PpSjlli = -iLUopVi° (72) 

PnrdrSflli + pprdrSfi°i = ~iuJopW° (73) 

6(3? = -liooil - e)v? (74) 

rdrS^f = -iujQ{l - e)w'} (75) 

Meanwhile the continuity equations lead to 

iuJor^Sp^ + drirpW^) - l{l + l)pVi^ = (76) 

and 

iuQpr^Sx? + dr[xp{l - Xp)rpw?] - Xp(l - Xp)l{l + l)/9u° + pWirdrXp = (77) 

Before we proceed, we will simplify the problem. Our aim is to determine analytic approximations for the funda- 
mental modes of the system, including the mutual friction damping. Solving the problem numerically is, of course, 
straightforward but does not lead to the same level of insight into the dependence on the various parameters. To 
facilitate an analytic solution, we will combine an incompressible background model with compressible perturbations. 
This simplifies the calculations considerably. In addition, since this is the same model that was considered by Lind- 
blom and Mendell [l3| we can compare our final results directly to the available literature. We thus assume that pn 
and p-p are both constant, while and 5pp are not. 

It is also useful to introduce a new variable for the co-moving degree of freedom. Let us define 

5hi ^ ^Spi = ^{pnSpi + PpSpl) (78) 

For a single barotropic fluid, dhi corresponds to the perturbed enthalpy. For a compressible background model we 
would have 



PnrdrSpl^ + pprdr-Spp = prdrShi - pdfiirdrXp 



(79) 
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However, for the uniform density model the gradient of the proton fraction vanishes so we simply have 

PnrdrSji'^ + pprdrSjip = prdrShi (80) 

It is also worth noting that 5hi has the same dimension as (5/3;. 
We now find that ([7^ and ([75)1 can be written 

Sh'^ = -iuJoVi° (81) 

and 

rdrSh^ = -iujaWl^ (82) 

Before we proceed, we need to decide what variables we want to work with. We can either remove [6hi,S/3i] or 
[(5/5;, (5xp] (or some other combination of these variables) from the problem using thermodynamic identities. Opting 
for the latter possibility, we use 

5pi = ^5hi + ^Spi (83) 

and 

54 = ^Shi + ^5f3i (84) 
In these relations we have defined, first of all, the speed of sound as 



We have also introduced 



and 



and made use of the identity [36 



dxp \ ^ 2 ( ^\ ^ 
dh)^ ^ { dp J I, {dp 

This reduces the number of required "thermodynamic" quantities to three; c^, ai and a2- 

For later convenience, it is useful to pause and consider the relative magnitude of the thermodynamic derivatives. 
To do this, take as an example an overall n ~ I polytrope with a proton fraction that is linear in the total density. 
This simple model is not completely unrealistic, and moreover it is straightforward to work out all the quantities we 
need. Assuming that 

p = Kp^ (89) 

we find that 

dp\ 2 { 



. , 2Kp = ci^[— =^ (90) 
Combine this with the assumption that the proton fraction (in equilibrium) is linear in the density. That is, take 

^) (91) 

PcJ 
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where a 10 ^ and pc is the central density of the star. This leads to 



dp\ _ 1 / dh 



and 



dl3\ 2Kp 

dxr- 



''P / h -^P -^P 

These estimates suggest that ai ^ a2 ^ Xp. Since we expect to have Xp ^ 1 it should be the case that 



(92) 



(93) 



(94) 



This agrees with the more realistic equation of state considered by Lindblom and Mcndell JJ| . We will make explicit 
use of this ordering later. 

Returning to the coupled system of equations, and combining the various relations we easily arrive at the two 
differential equations 



and 



drir^drSh^) -1(1 + 1) 



drir^drSp^) - 1(1 + 1) 



1 - 



/(/+I)c2 

(1 - e)u!Qa2r'^ 



(95) 



Z(/ + l)(l-.Tp)c2 



6(3\ 



, (l-£Voai^%,o 



3^p(l 2;p)c;, 



= 



(96) 



These are the equations to be solved. 

Before we proceed, it is useful to introduce dimensionless variables. First we introduce loq = aau) where lu'^ 
GM/R^. Then we consider a new radial variable 



(97) 



This means that have 



ds{s^ds5h!l)-l{l + l) 



1(1 + 1) 



6h°i + ai(jls^50{ = 



together with 



ds(s^ds5(3f)-l(l + \) 



1 - 



(1 - e)a2<7ls'^ 



1(1 + 1) 



(1 — £)q;i(TqS^ 



Sh'} = 



(99) 



For simplicity, we have assumed that Xp is small (^ 1). This should always be the case. When the equations are 
written in this form it becomes apparent that the coupling term in ([55]) is more important than that in (|98p . From 
this one can deduce that there should exist solutions to the problem such that Sh'i ^ SPf. These are the modes that 
we will focus on. This is natural, if our main focus is on oscillations that radiate gravitational waves at a significant 
level, e.g. by being driven unstable 0]. 

As already mentioned, since we have assumed that the background configuration is uniform, our model is identical 
to the incompressible/compressible model considered by Lindblom and Mcndell [T^ . From their work we know that 
we can write down the solution to the coupled equations in closed form using (spherical) Bessel functions. This 
solution would contain all the modes of the system, fundamental modes and pressure modes with varying degree of 
CO- and counter-moving character. However, this solution is not very practical for our present purposes. If we want 
to solve the order Q problem explicitly, rather than estimate the mutual friction damping via the energy integral 
approach (as Lindblom and Mendell did), we need to be able to solve another system of equations where the leading 
order mode-solution acts as source. Expressed in terms of the Bessel-function solutions, the order problem is very 
messy. Hence, we opt for a different strategy and introduce yet another simplifying approximation. 
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In order to proceed analytically, let us assume that <C 1, in which case we can attempt to solve the problem 
using a power series. Is this reasonable? Well, let us again consider the case of an n = 1 polytrope. In that case 



K = 2nG i^-j (100) 
and it follows that 

where p is the average density of the star. This shows that, in our uniform parameter model, the power series Ansatz 
makes sense as long as we assume p 3> 3p. That is, the calculation should be relevant for the conditions in a neutron 
star core. However, it is obviously not completely consistent. The assumptions will not hold near the surface of the 
star, since one tends to have — > as r — i?. However, since the surface region is already dealt with in a rather ad 
hoc way this does not concern us too much. 

Now that we have a small parameter in the problem, we can try to find a power series solution. It is natural to first 
rewrite the coupled problem as a single fourth-order equation for (say) 6h'^. This is easily done by combining (j98p 
and Making the Ansatz ^ 

N 



6h^ ^s^J2 (102) 



we find that the first few coefficients are determined by 

ai = — az ~ (103) 

and 

-^^^ ^-^ao-2(l + a2)(2/ + 3)a2--^ ^ -a4 = (104) 

Xp (Tq 

We now insert this solution into (l98l). If we write 



N 



<5/3? ^s'J2 bn-s"" (105) 



n=0 



then we must have 



and 



,^_<^iao±2(2^ (106) 



^ _^tl±^^±3^ (107) 
To complete the solution, we need to satisfy the boundary conditions. We want 

djpf = at s = — (108) 

Cs 

Keeping the first two terms in the series for Spf, this condition leads to another relation between the three coefficients 
ao, a2 and 04. Combining this relation with p04p we arrive at an expansion for Sh'j' where the overall scaling is given 
by ao, and the only other unknown parameter is the frequency ctq- To fix the frequency, we impose the remaining 
boundary condition. That is, we require 

drSh'^ + Sh\' = at r = R (109) 
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where (for a uniform background model) 



AnGpR 



(110) 



This leads to the condition 



LOJx 

Some algebra now leads to a solution with frequency 



That is, we have 



21 



at 



QR 



(111) 



(112) 



IGM 



21 



LUR 



(113) 



The leading order result is exactly what one would expect for an incompressible fluid ball in the Cowling approximation. 
The first correction to this is directly associated with the compressibility. The presence of the second fluid degree 
of freedom, e.g. the link to (5/3, appears at the next order of approximation. It is also worth pointing out that the 
solution is such that (omitting geometric factors of I) 



bo 



LORV 
Cs ) 



(114) 



This demonstrates that the mode we have determined is such that 5h^ ^ 5(3^ . The associated fluid motion is, indeed, 
predominantly co-moving. The conclusion that the co-moving f-mode is hardly at all affected by the two-fluid nature 
of the system accords well with the results of Lindblom and Mendell [l^l ■ 



B. The slow-rotation corrections 



Having approximated the f-mode solution to the non-rotating problem, we will now work out the first order slow- 
rotation corrections. This will include the mutual friction damping. 
The equations that need to be solved at order are, first of all 

2{l + 2){W^ - ZF,")Q,+iy,™ - 2{l - 1)[< + + l)V^']QlY,^^ + l{l + l)oJoUlYr = (115) 

which determines the axial rotational correction Ul to the f-mode. We are not going to determine this quantity here 
because it does not affect the mode damping, which is our main concern. A similar equation for the counter-moving 
degree of freedom determines the axial correction, u\. This is also not of immediate relevance for our discussion, so 
we do not consider it. 

From Section IV we see that the equations we actually need to solve are 

+ l)cji - 2m]Vi^ + l{l + l)uJoVi^ - 2mW^ - + l)5h] = (116) 

and 

irdrSh] + 2mV[' - uJiW^ - loqWI = (117) 

together with 

- + 1)513] + 2mB'w^ - [l{l + l)tJi(l - e) - 2mB']v^ - l{l + l)tJo(l - e)v] 

+ 2iB[\ -{1 + 3)Qf+, + 2)Qf]w'^ + 2iB[m^ + l{l + 3)0^+1 + {I + - 2)Q2]„o|y™ 

- 2i{l + 3)(lt;° - w°)eg,+iQ/+2i"z+2 - 2i(; - 2)K - + l)v°i]BQiQi-iYi"^^ = (118) 



15 



and 

[irdr6Pl - [oj.il - e) - 2iB{\ - Qf ~ Qf+,]w'i' - u;o{l - e)w] + 2[mB' - + 1)0^ ^ ^Qf^ JB]i;°}rr 

+ 2i(h>? - w'})BQi+iQi+2Yi"l^ - 2i[w°i + {I + l)v1]BQiQi^iYi^22 = (119) 

If we want to determine the rotational correction to the f-mode, then we only need to consider the order 11 terms 
that are sourced by non-rotating terms. Thus the problem reduces to solving 



and 



(120) 
(121) 



irdrSh] - uj^Wl = uJiW^ - 2mVl^ 
together with the continuity equation (assuming a uniform background) 

iuJar^Spj + pdr{rWl) - l{l + l)pVi^ = -iujir^Sp^ (122) 
Note that there is no multipole coupling in these equations. For the other degree of freedom we have 

l{l + l)[iS(3l - ujq{1 - e)vl] = -2mB'w'j' + [l{l + l)tJi(l - e) - 2mB']v^ 

- 2iB[l + 3)gf+i + (/ - 2)Qf]w? - 2iB[m^ + l{l + 3)Qf+i + {I + - 2)gf]w;" (123) 

and 

irdrSfil - cjo(l - e)w] = [wi(l - e) - 2iB{l - Qf - Q'i+{)]w°i - 2{mB' - + l)Qf - lQf+i]B}v^ (124) 
together with 

iujor^dxj + .Tp(l - Xp)[drirwl) - l{l + l)v]] = -iwir^fe^* (125) 
We want to find solutions that satisfy the boundary conditions 



at r = i? 



(126) 



and 



iujodrSPl = ~iujidr50{ = , at r = i? (127) 
After some manipulations (making use of the leading order relations) we arrive at the two coupled equations 



dr(,r'^dr6h\)-l{l + l) 

and 

dr{r^dr5l3l)-l{l + l) 



,2J2 



1 



/(/+I)c2 

(1 - £)ujla2r'^ 



Shi + d^spl = -^^^ 



(128) 



;(Z + l)(l-a:p)c2 
2woQ!ir-^ 



,2^2 



6^1 



3^p(l Xp)Cg 



Xp{l 2;p)c^ 
2iB 



[{1 - e)uji - iB{l 



(1 - e)ujo 



[l + {2l-l)Qf-i2l + 3)Qf+,]rdrSp^ 



2iB 



(1 - e)ujo 



[m' - l{l + + + l){2l - l)Qf + 1(21 + 3)Qf+,]Sp° 



2a2Wo'"^ 
(1 - .Tp)q 



■[(1 - £>i - t{l -Qf- Qf+,)B]S/3f (129) 
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The problem has the anticipated form, a coupled system of equations for Shj and SPl which differs from the non- 
rotating problem only by the presence of leading order source terms. To proceed, we will follow the same strategy as 
in the leading order calculation. However, it is beneficial to first note that the source term in the second equation 
simplifies considerably if we focus on the I = m f-modes. Since Q^+i = l/(2m + 3) we see that the factors in front 
of the first two SPi pieces in the right-hand side of (|129[) are then identically zero. Moreover, we know from the 
leading order calculation that 5/3° <^ 6h^, cf. ()114p . As long as we are only interested in the leading order rotational 
correction and the leading order mode damping, this allows us to neglect also the remaining 6(3i part of the source. 
Using ujqt/cs = <Jos as well as ui = uji/ujq and assuming that Xp ^ 1 we then arrive at the simplified equations 



ds{s'ds5h\)-l{l + l) 



and 



t2 o2 



dsis'dsSpD-iii + i) 



1 - 



(1 - e)(Tla2s'^ 



1(1 + 1) 



8h\ + a^ais^SPl = -2alais'^Sh^ 



(130) 



(131) 



where 



(1 - s)ai 



2i{l + l)B 



aQU){2l + 3)_ 

Combining the two equations, making the same Ansatz as in the non-rotating case 



(132) 



N 



5h] = s' Qni 



n=0 



and taking the source term to be Sh'^ = Cs^ |4T| , we find that 



ai = as 



(133) 



(134) 



Meanwhile, we have 



a-'^ial - a2Xp) 



aa~2{l + a2){2l + 3)a2- 



8(2Z + 3)(2; + 5) 



-04 



- ao<Tia2 C = 



(135) 



The difference now is that we are only interested in the particular solution due to the presence of the source term. 
In order to remove the unwanted homogeneous solution we set ao = 0. Then ()135p becomes a relation between the 
known mode amplitude C and the two coefficients 02 and 04. As in the non-rotating problem, we get a second such 
relation from the boundary condition for S^f. From (jlSOp and 



Sf3l 



N 
n=0 



(136) 



we get 



ai 



(2^ + 3) ■ 
2 ^2 + CTiG 



(137) 



and 



0,2 



4(2/ + 5)04 



(138) 



This solution has to satisfy the surface condition 



dsSPl =0 at s = — 



(139) 
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This leads to another relation between C, 02 and 04. This means that we can write down the solution for Shj with an 
overall scale C (as expected for the particular solution) and with the frequency correction cti as the only remaining 
undetermined quantity. 

The final condition to be satisfied can be written (again, for / ~ m) 



ds5h\ - ^ol5h] 



u)R 



ai [al + I) 



21 
ujao 



Cs' 



at 



(140) 



Inserting the power series solution, one can show that this condition leads to the leading order rotation correction 
of the f-mode frequency being 



Re cTi 



1 



Meanwhile the leading damping term (the imaginary part of cti) is 



Im (Ti 



Vl{2l + 3)^(21 + 5) Xp 



(141) 



(142) 



These are the final results of the f-mode analysis. After retracing our steps to recall the various definitions, we find 
that the mutual friction damping follows from 



Im uji 



2{l + l){3l + 5) al [GmV 2(/ + l)(3/ + 5) 1 fdpV [GAiy 

{2l + 3r{2l + 5)xl{Rcj) i2l + 3n2l + 5)pl\d(3)^{ R J ^ ' 



That is, the damping timescale is 



1 



_ (2^ + 3)3(2^ + 5) 
Im ujiQ ^ 2(Z + l)(3Z + 5) 



1 fdpV fGM"- ^ 



pI R 



1 



(144) 



This result completes our analysis of the (co-moving) f-mode in a supcrfluid neutron star. 

It is obviously relevant to compare the estimated f-mode damping timescale to previous work. To do this, we first 
need to recall that all previous work has focused on the case where electron scattering off the vortex array is the main 
cause of mutual friction. Then we have [13, uM 



1/2 



0.05 



7/6 



1/6 



10i4g/cm'' 



(145) 



where we have used the relation between the entrainment and the effective proton mass; 



1- ^ 



(146) 



Taking rrip/mp ~ 0.3 we have B « 5.5 x lO^'* in good agreement with the result used by Lindblom and Mendell 
in their investigation of the f-mode problem [iJl- The overall scaling with I in (|144p also appears to be similar to 
their result. This is evident from the results in Table [T] which compares our results to data from Table 1 in p^ . 
This comparison shows that, in the m*/mp = 0.3 case, our damping times are about a factor of 2 longer than those 
estimated by Lindblom and Mendell. The results also differ in the predicted dependence on the entrainment. In our 
calculation, the entrainment only enters (|144p indirectly through its effect on B. The data given by Lindblom and 
Mendell hints at a different behaviour. In the m*/mp ~ 0.8 case we find that the difference between our damping 
result and Table 1 in [l3| is closer to a factor of 3. 

Most likely the main difference originates from the use of the energy integral approach in one case and the direct 
determination of dissipative mode solutions in the other. In order to check this assumption, we have estimated the 
mutual friction damping using the energy integral approach (together with our leading order f-mode solution). That 
is, we evaluate p6p for the non-rotating f-mode solution. Then we find that 



(2^ + 3)3(2/ + 5) 
6(2P + 6^ + 5) 



1 (GMV 



d(3j^y R J 



1 

Bn 



(147) 
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TABLE I: Estimated mutual friction damping timescales for the (co-moving) f-mode of a superfluid neutron star. We compare 
our results to previous work by Lindblom and Mendell [TJ]. The parameters are those discussed in the main text, and correspond 
to nip/mp = 0.3. The stellar model has radius 15 km and average density 4x 10^'' g/cm"^. Note that, in our case the non-rotating 
f-mode frequency is given by = 4ZnQ/3, where Qo = y/TvGp and p is the average density. The difference between the listed 
frequencies (15% for 1 = 2 and decreasing with increasing I) should mainly be due to our use of the Cowling approximation. 
The dissipative mode result (|144|l for the mutual friction timescales (which scales as Qq /Q) differs by about a factor of 2 from 
the tabulated Lindblom-Mendell results. By comparing to (|147p we learn that this is mainly due to different geometric factors 
in the energy integral approach and the direct dissipative mode calculation. 
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This scaling with I here differs somewhat from that in p44p . This introduces a numerical factor of w 2 for small values 
of I. This factor brings our results very close to those of Lindblom and Mendell, cf. Table [D The main difference 
between the energy integral result and the full dissipative mode calculation is a geometric factor of order unity. This 
is what one would expect. At the end of the day the astrophysical implications of the results are the same. 

Before we move on, it is worth making a comment on the apparent lack of entrainment scaling. The damping 
timescale in (|144p does not depend directly on the entrainment parameter e due to a scries of cancellations that 
occur when we impose the order D, boundary condition Vidf3 = 0. If we were to impose the (slightly more realistic) 
condition 5w^ = at the surface (corresponding to a common surface), these cancellations would not occur and there 
would be an explicit dependence on e in the damping timescale. This does not, however, affect the numerical results 
significantly. 

C. The f-mode instability window 

To conclude the discussion of the superfluid f-mode, it is worth considering the impact of the results on the 
gravitational-wave driven instability of this mode. The close agreement between our mutual friction damping rates 
and the results of Lindblom and Mendell [l3| obviously means that their key conclusions stand. That is, the instability 
of the f-mode is likely to be completely suppressed in a superfluid neutron star. However, we think that this result has 
sometimes been misunderstood. The result does not show that the secular f-mode instability cannot play a role for 
astrophysical neutron stars. To show this, we have combined the different timescales for gravitational-wave growth of 
an unstable mode with the damping due to shear- and bulk viscosity from Ipser and Lindblom [3^ . The results, for 
the I = m = 4 f-mode that leads to the strongest instability in a Newtonian star, are shown in Figure [TJ The data 
in the figure corresponds to an n = 1 polytrope with mass 1.5Mq and 12.533 km (the average density is 3.6 x 10^'* 
g/c m'^). In order to connect this with our mutual friction approximation, we have used the model parameters from 
i.e., 

4 X 10^'*g/cm^ 
0.06 

1.911 X 10"VVcm'^ 

Combined with the canonical value for B given in Ijl45p these parameter values lead to the mutual friction damping 
completely overwhelming any gravitational- wave driving of the f-mode. Table [T] provides similar results for the model 
considered in [l^ . a star with radius 15 km and average density 4 x 10^^ g/cm^ (which means that the mass is quite 
large, 2.84Mq). We have scaled the frequencies and the timescales using Q,o — s/i^Gp where p is the average density. 

However, this result is only relevant below the critical temperature at which the stars core become superfluid. 
Suppose that we take the critical temperature to be 5 x 10^ K, which is a typical value [l^. Then the f-mode 
instability window remains unaltered in hotter stars. This is evident from Figure [TJ Of course, as soon as a sizeable 
part of the core is superfluid, the instability is no longer present. However, it seems that there is still scope for the 
unstable f-modes to play a role in the evolution of a nascent neutron star born spinning near the breakup velocity. One 
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FIG. 1: The f-mode instability window for the I = m = 4 f-mode. The data for gravitational radiation reaction, shear- and bulk 
viscosity are taken from [ssj ). The instability is active above a critical rotation rate (thick solid line) at each given temperature. 
The mutual friction, which acts only below the superfluid transition temperature (here taken to be 5 x 10® K, lower temperatures 
are indicated by the grey region in the figure) is estimated using (|144p and the parameter values given in the main text. To 
illustrate the role of a weak mutual friction we show (as thin solid lines) the instability curves for B in the range 10^^'^ — 10^^, 
more than three orders of magnitude weaker than the canonical value (|145|l . The solid rotation rate of the star is given as a 
fraction of the breakup rate Qk w 0.639f2o . 



should also remember that the instability is stronger in a relativistic model. In fact, in this case the / = to = 2 mode 
may also become unstable. Based on the available evidence it vifould be premature to rule out the f-mode instability 
for realistic neutron star models. The problem requires further attention. 

In Figure [T] we also show the effect of a weaker mutual friction. Suppose that (|145p is. for some reason, not the 
typical value. Our understanding of neutron star core physics is not complete, so it is interesting to consider a range 
of possibilities. The different thin solid curves in the figure show the effect of mutual friction for the given values of 
B. For the considered 1.5Mq model the results correspond to (for / = 4) 

Tilo « 2.5 X 10^ (^) (-^) (148) 



From these results we learn that the mutual friction must be at least three orders of magnitude weaker than the 
canonical value in order for the f-mode to be unstable below the superfluid transition temperature. 

Finally, it is worth noting the following. The coefficient B' was not present in our final f-mode equations. This 
is due to a series of, perhaps surprising, cancelations. The upshot of this is that our results arc also valid in the 
strong drag regime. Expressed in terms of the drag parameter TZ, the data in Figure [1] shows that the mutual friction 
suppresses the f-mode instability in a superfluid neutron star with 10"^ < TZ < 10^. This conclusion is interesting 
since the strong drag regime has not been considered before. It also shows that the suppression of the f-mode takes 
place for much of the plausible parameter range. 



VII. THE R-MODES 



Having explored the f-mode in a superfluid star, revisiting the issue of the mutual friction damping, we will now 
consider the Coriolis driven r- modes in a slowly rotating star. The r- modes are interesting because they also suffer 
a gravitational- wave driven instability In contrast to the f- modes, which only become unstable at fast rotation 
rates, the r-mode instability may set in already at quite modest spins. 
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In a single fluid star, the r-modes are purely axial to leading order. Moreover, their frequency is linear in the 
rotation rate. Hence, it is natural to make the Ansatz 

uj = uj^n (149) 

together with 

Wi = VLWl , Vi = nVi^ , Ui = Uf (150) 

and 

wi — Qwl , vi ~ rivl , ui ~ riu^ (151) 

Note that, if we want to work out the order fl corrections to the mode (with the above ordering) we will first need 
to account for the centrifugal force and the change in shape of the star. We will discuss this problem elsewhere [1^. 
Here we will focus on the problem at linear slow-rotation order. 

From the general slow-rotation equations in Section 4, we immediately see that the r-mode assumption decouples 
the two degrees of freedom. First of all, the average vorticity equation leads to 

+ 1)^0 - 2m] C/f Y," = (152) 

This shows that we must have a single multipole solution, with frequency 

2m 

To determine the associated eigenfunction we consider the divergence equation and the radial Euler equation. These 
lead to the recurrence relations 

- in{n + l)dhl - 2(n - l)(n + l)QnU°_, - 2[{n + i f - l]Qn+iU°+, = (154) 

and 

irdrShl + 2{n - l)QnU^-i - 2(n + 2)Q^+iU°+^ = (155) 

For simplicity, we have assumed that the background is uniform (as in the f-mode analysis in the previous section). 
These equations show that the only way to have a single multipole axial solution is to have ^ 0. This follows since 
Qm = 0. In other words, we will have non-trivial modes only for I = m. Inserting n = m — 1 in the two equations we 
have 



and 



These combine to 



and the familiar solution 



t(m + l)Shl^+, - 2mQ,„+i[/0 = (156) 

irdrShl^+-^ + 2mQ,„,+iU^:^0 (157) 

rdrU°, - {m + 1)U°„ = (158) 

C/O = Ar"+i (159) 



This analysis shows that, to leading order, the standard r-mode remains unchanged in a superfluid star. We also see 
that we need to go to higher orders in rotation if we want to determine the mutual friction damping of these modes. 
Such calculations have been carried out by Lindblom and Mendell [1^1 and Lee and Yoshida [l^|. Motivated by the 
recent evidence that the strong drag regime may be relevant, we are currently revisiting this problem [36j . 
Now consider the counter-moving degree of freedom. In that case, the difference vorticity equation leads to 

{l{l + 1)(1 - e)uJo - 2mB' - 2i[l{l + 1) - m^]B} u^Y;" = (160) 
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That is, there should exist a single multipole solution with frequency 

As in the co-moving problem, the associated eigenfunctions follow from the divergence equation and the radial Euler 
equations. These lead to the recurrence relations; 

- in{n + 1)5/3^ ~ 2(n + l)[(n - 1)B' + mi6]Q„<_i - 2n[{n + 2)B' - imB]Qn+iu'i+i - (162) 

and 

irdrSpl + 2[{n - 1)B' + im^]g„<_i - 2[{n + 2)8' - imB]Qn+iul+i = (163) 

Again, it is easy to see that the only way to have a single multipole solution is to have 7^ 0, i.e. we must have 
I = m. This leads to 

- iim + l)5/3i+i - 2m[B' + iB]Q„,+iu"^ = (164) 

and 

trdrSPl+i + 2m[B' + i6]Q„+i< ^ (165) 

That is, we have 

rdrul^-{m + l)ul^^Q (166) 
which means that the counter-moving solution also takes the form 

= Br™+i (167) 

What do we learn from this exercise? First of all, we should recognize that we have been somewhat cavalier in 
the discussion. Since we have assumed that e, B' and B are all constant, the analysis leading to (|16ip is clearly only 
valid for uniform background models. This tells us that the purely axial counter-moving solution only exists for this 
simplified model. In a more general case, this mode will become an axial-led inertial mode. In the weak drag regime, 
these inertial modes have been determined numerically by Lee and Yoshida [l^. 

The counter-moving r-modes are nevertheless interesting. Two particular features are worth noting. Let us first 
consider the mode pattern speed 

<7p 168 

m 

In the case of the normal r-mode we see from (|153p that the pattern speed is always negative. That is, these modes are 
retrograde with respect to the star's rotation. This is the criterion that renders the mode unstable to gravitational- 
wave emission at all rotation rates (in an otherwise non-dissipative star). Meanwhile, from (|161[) we find that 



2m -, 1 2m 

-B' = — 



1 — e (m -I- 1) 1 — En/xp (m + 1) 



1-^ 



(169) 



This relation shows that the second class of superfluid r-modes may, in fact, be prograde. For this to be the case we 
must have (assuming that 1 — e > , see below) 

B' > Xp (170) 

which may well happen. Recall that the mutual friction parameter is related to the induced friction on the vortex. 
From (HI we see that we need to have 

> — - — > 1 (171) 

1- Xp 

in order for the mode to be prograde. Clearly, systems in the strong drag regime (where TZ — > 00) will satisfy this 
condition. Alternatively, we may require 

— >1 (172) 
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Recalling that £„ = (pp/pn)ep and using the relation between the entrainment and the effective proton mass (|146p . 
we find that (|172p corresponds to 

in* 

Xp>^ (173) 
nip 

This condition is unlikely to be satisfied in a neutron star core, where typical values would be Xp « 0.1 and m* « 0.5mp, 
but the possibility is nevertheless interesting. In particular since the mode will actually be unstable (due to the presence 
of the mutual friction) if the condition is met. Since S > it is clear from (jl6ip that the imaginary part of the mode 
frequency is negative if (|172p is satisfied. The existence of this unstable regime is interesting, at least conceptually. 

As a final check let us compare the damping timescale calculated in (|16ip with that calculated using the integral 
approach, cf. p6|) . Using the definitions for the energy, p4)) . and the dissipation due to mutual friction, ((33)) . one 
readily finds for the counter-moving r-mode solution; 



dtEe = -2p{l - Xp)Bn / [l{l + 1) - m'^]\u'^\^r'^dr 
Jo 



(174) 



Meanwhile 



These lead to 



1 



R 



E ^ -pxpil - xp){l - e) I lil + l)\u^\'r^dr (175) 







Imcj = - = - r , (176) 

r (l-£-) 1(1 + 1) ^ ' 

which agrees perfectly with the damping timescale extracted from p6ip . Since we are using the full dissipative 
mode-solution in the energy integrals, this is as expected. 



VIII. CONCLUDING REMARKS 



The aim of this paper was to lay the foundation for a renewed assault on the problem of dissipative superfluid 
neutron star oscillations. We have discussed the oscillations of slowly rotating superfluid stars, taking into account 
the mutual friction force at linear order in the (presumed) slow rotation of the star. We have considered both the 
fundamental f-modes and the inertial r-modes. Our analysis differs from previous studies in that we do not assume 
weak mutual friction from the outset, the final results are also valid in the strong drag regime. 

In the case of the f-modes, we worked out an analytic approximation for the mode which allowed us to write down 
a closed expression for the mutual friction damping timescale. This result, which is in good agreement with previous 
numerical results of Lindblom and Mendell [l4|, provides useful insight into the dependence on, and relevance of, 
various equation of state parameters. The scaling with the harmonic index I is also obvious from our final formula. 
The analysis is readily extended to stars with superfluid cores that do not extend all the way to the surface (as 
assumed in our analysis), although the result is then less transparent. 

In the case of the r-modes, we have confirmed the existence of two classes of modes. However, we demonstrated 
that only one of these sets will remain purely axial in more realistic situations. This agrees with previous results of 
Lee and Yoshida [l^. We discussed some peculiarities of the counter-moving r-modes. In particular, the fact that 
they may be unstable for some parameter values. Even though we do not expect this instability to be relevant for 
realistic superfluid stars, its existence is of conceptual interest. 

Building on the formalism and the results presented in this paper, we are currently carrying out a detailed study 
of the mutual friction damping of the r-modes at second order in the slow-rotation approximation [3^ . At the same 
time we are considering neutron stars with exotic hyperon and/or quark cores. Since the multifluid aspects of those 
problems have never been considered in detail, these are exciting developments. They are, in fact, necessary if we 
want to understand the dynamics of realistic models of mature neutron stars. 
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